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Abstract 


This paper has developed a complete two-phase model of a proton exchange membrane (PEM) fuel cell by considering fluid flow, heat 
transfer and current simultaneously. In fluid flow, two momentum equations governing separately the gaseous-mixture velocity (u,) and the 
liquid-water velocity (uw) illustrate the behaviors of the two-phase flow in a porous electrode. Correlations for the capillary pressure and the 
saturation level connect the above two-fluid transports. In heat transfer, a local thermal non-equilibrium (LTNE) model accounting for intrinsic heat 
transfer between the reactant fluids and the solid matrices depicts the interactions between the reactant-fluid temperature (Tp) and the solid-matrix 
temperature (7,). The irreversibility heating due to electrochemical reactions, Joule heating arising from Ohmic resistance, and latent heat of 
water condensation/evaporation are considered in the present non-isothermal model. In current, Ohm’s law is applied to yield the conservations in 
ionic current (im) and electronic current (i,) in the catalyst layer. The Butler—Volmer correlation describes the relation of the potential difference 


(overpotential) and the transfer current between the electrolyte (such as Nafion™) and the catalyst (such as Pt/C). 


© 2006 Elsevier B.V. All rights reserved. 
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1. Introduction 


As illustrated in Fig. 1, three kinds of two-phase transport 
phenomena can be characterized in an operating proton 
exchange membrane (PEM) fuel cell. The first one is two-phase 
flow, i.e., Uw and ug. It is attributed to the morphological phase 
in water. Pressure difference between the liquid water and the 
gaseous mixture causes the capillary force that affects the satu- 
ration level in the porous electrode. The second one is two-phase 
temperature, i.e., Tg and T;. Temperature difference between the 
reactant fluids and the solid matrices (i.e., carbon fibers) results 
in intrinsic heat transfer in the gas diffusion layer. The last one 
is two-phase current, i.e., im and ig. Potential difference between 
the electrolyte (such as Nafion™) and the catalyst (such as 
Pt/C) drives the transfer current inside the catalyst layer. 

A large amount of research works about two-phase models of 
a PEM fuel cell have been proposed in the past few years. Only 
those most closely related to the present study are briefly cited 
below. Berning and Djilali [1] and Nam and Kaviany [2] mod- 
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eled two-phase flow in a PEM fuel cell. Two-phase currents were 
not considered since the catalyst layer was treated as a thin inter- 
face. Zhou and Liu [3], Ju and Wang [4], and Lin and Nguyen 
[5] coupled two-phase current to two-phase flow under the adi- 
abatic conditions. Hwang et al. [6-9] developed a local thermal 
non-equilibrium (LTNE) model to discuss the characteristics of 
two-phase temperature in a PEM fuel cell. The reactant fluids 
were assumed to be gaseous. Basically, all previous two-phase 
models discussed one or two of the aforementioned two-phase 
transports. None has provided a full coverage of the three two- 
phase transports so far. The objective of this work is to develop 
a complete two-phase that covers all two-phase components 
in a porous cathode of a PEM fuel cell, i.e., velocity (liq- 
uid water/gaseous mixture), temperature (reactant fluids/solid 
matrices), and current (ionic conductor/electronic conductor). 
This paper is the first step in developing a comprehensive two- 
phase model for a conventional PEM fuel cell cathode, which is 
believed to be of great importance in accurately describing the 
transport phenomena occurring in a PEM fuel cell. 


2. Model description 


A schematic drawing showing a two-layer cathode of a PEM 
fuel cell is given in Fig. 1. The humidified air diffuses into the 
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Nomenclature 

c mole concentration of the species (mol m7?) 

D; diffusivity of species i (m? s—!) 

F Faraday’s constant (96,487 C mol`!) 

hy interfacial heat transfer coefficient (volumetric) 
(W m? KT!) 

i current density (A m~?) 

j transfer current density (A m°) 

k thermal conductivity (W K=! m7!) 


Mphase Phase change mass transfer (kg m~? s7!) 

MORR,O, OXygen mass transfer rate due to oxygen reduc- 
tion reaction (kg m~? s7!) 

MORR, H0 Water mass transfer rate due to oxygen reduc- 
tion reaction (kg m~? s7!) 


M molecular weight (kg mol~!) 
p pressure (Pa) 

d heat flux (W m7?) 

R universal gas constant (W mol`! K7!) 
s saturation 

T temperature (K) 

u velocity vectors (m s7 1) 

Xj mass fraction of the species j 
x,y coordinate system, Fig. 1 (m) 
Greek symbols 

a symmetric factor 

E porosity (gas diffusion layer) 
0) phase potential (V) 

K permeability (m7) 

p density (kg m?) 

o conductivity (Q7! m7!) 

T tortuosity 

w mass fraction 

Subscripts 

c capillary or catalyst layer 

eff effective 

f fluid 

g gas phase 

int intrinsic 

m membrane (electrolyte) phase 
ORR oxygen reduction reaction 

s solid or catalyst phase 

w water 

Q resistive loss 

Superscript 

sat saturation 


gas diffusion layer (GDL) and catalyst layer (CL) from the flow 
channel. Then, the oxygen reacts together with the protons to 
form liquid water in the catalyst layer. Meanwhile, the heat gen- 
eration due to the overpotential and irreversibility heating during 
the electrochemical reaction is removed by either the reactant 
fluids or the solid matrices. The porous cathode is treated as 


homogeneous porous media with uniform morphological prop- 
erties such as porosity and permeability. Typical data are list in 
Table 1. 


2.1. Two-phase velocity 


Two-phase flows in porous electrode follow separate equa- 
tions for the gaseous mixture and liquid water: 


V(PgUg) = MORR, O, — M phase (1) 


V(OwUw) = MORR,H2O + Mphase (2) 


The first term on the right-hand side of above equations means a 
sink/source of oxygen/water during the oxygen reduction reac- 
tion (ORR). Since the reacting temperature of oxygen reduction 
reaction in the cathodic catalyst layer should be below 100°C, 
the product water (79RR,H,0) is assumed to be in liquid phase. 
Through non-equilibrium evaporation (—/i1phase) of liquid water 
from the pre-existing liquid surface, water vapor is introduced 
diffusively into the porous medium. pg and pw are the densities 
of gaseous mixture and liquid water, respectively. The veloci- 
ties of gaseous mixture and liquid water are expressed by the 
Darcy’s law, i.e., 


£[LUg 


1 
-V - (pUglg) = —EV Pg + V - (Ug VUg) — (3) 
E Kkrg 
1 
ae -(OUyUw) = —ELV Pc + (Pw — Pg)gl] + V - (UwVUw) 
= Eww (4) 
Kkrw 


where xo is the permeability of the dry electrode, while krg and 
kw are the relative permeabilities of the gaseous mixture and the 
liquid water, respectively. They account for the reduction in pore 
space available for one phase due to the existence of the second 
phase. In this work, we assume isotropic porous electrode with 
relative permeability of each individual phase preoperational to 
the cube of phase saturations, i.e., krg = (1 — sy and kw = s° [10]. 

The capillary forces due to the imbalance of the liquid-water 
pressure and the gaseous-mixture pressure can be expressed as 


Pc = Pw — Pg (5) 
This shows the liquid-water pressure should be higher and makes 
it possible to transport liquid water out of the hydrophobic 
porous electrode without severe flooding, and also prevents 
wicking of water from the channels. The functional variation 
of capillary pressure with saturation is prescribed following 
Leverett [11]: 


E 
Pc(S) = Pw — Pg = ô (=) f(s) (6) 
Ko 
where ô is the interfacial liquid/gas tension; € the porosity; f(s) 
is the empirical function proposed by Udell [12]: 
f(s) = 1.4201 — s) — 2.121 — 5)? + 1.2601 — s)? (7) 


The interfacial mass-transfer rate between the gas and liquid 
phases 7itphase is assumed to be proportional to the amount of 
the reactant in the porous media and the difference between 
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| Components | Gas diffusion layer (GDL) 


| Catalyst layer (CL) | 


| Velocity | Two-phase, Ug, Uy 


Two-phases, Ug, Uy 


Temperature | Two-phase, Ts Tr 


| Single-phase, TT r | 


Current Single-phase, i, 


y, mm 
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0.35 


0 Flow Channel 0.5 


Two-phase, is, im 


Rib Shoulder 1.0 


Fig. 1. Schematic drawing of model domain of a two-layer porous cathode. 


water-vapor partial pressure (py,0) and its saturation pressure 
(pity), i.e. [2,13], 


(Pmo— Pio) 


kcone(1—s)xH,0 RT; 


í t 

if pmo— Piho Z 0 

Mphase = p 
Ww sat 

kevpEs Mao (Pmo— Pho) 


i t 
> if pmo — Pio <0 
2 


(8) 


The upper half represents the condensation rate while the lower 
half represents the evaporation rate. kcon and keyp are the con- 


densation and evaporation rate constants, respectively [13]. xH,0 
is the mole fraction of water vapor. The saturation pressure of 
water, pio is a function of the operating temperature (Tẹ), 
which is read directly from the steam tables. 


2.2. Two-phase temperature 


A local thermal non-equilibrium model [8] is used to describe 
the two-phase temperatures in the porous electrode. The energy 
equations for solid matrices and reactant fluids in the GDL are, 


Table 1 

Porous-electrochemical properties of the present modeled fuel cell 

Name Description Value 
COp ref Reference mole fraction of oxygen (mol m~?) 3.6641 
CH20,ref Reference mole fraction of water (mol m~?) 0.0703 
ig Exchange current density (A m7?) 5.0 x 10-4 
Sa Surface area to volume ratio (m~!) 1000 
Jo (Volumetric) exchange current density (A m~?) 0.5 

E Porosity 0.45 

Ko Absolute permeability (m?) 2.55 x 10-18 
hy Intrinsic heat transfer coefficient (W m~? K7!) 1x 10° 
Om Tonic conductivity (electrolyte phase) (Q7! m7!) 14.4 

Os Electronic conductivity (catalyst phase) (27! m7!) 300 

Kes Solid-phase thermal conductivity in the catalyst layer (W kg~! K7!) 1.1 

Kg eff Solid thermal conductivity (W kg! K~!) 17 
Keett Fluid thermal conductivity (W kg! K`!) 0.051 
ke,eff Catalyst layer thermal conductivity (W kg™! K7!) 0.656 
kı Electrolyte thermal conductivity (W kg—! K7!) 0.5 

Cpt Reactant-fluids thermal capacity (Ikg~! K7!) 1006 
Cp,s Solid thermal capacity (JIkg~! K7!) 885 

Ps Solid-part density (kg m~>) 1100 
pi Electrolyte (Nafion 117) density (kg m~’) 1960 
Pw Water density (kg m~?) 1000 
lly Water viscosity (Pa s) 3.65 x 1074 
EGDL Porosity of gas diffusion layer 0.48 
ECL Porosity of catalyst layer 0.42 
kevp Evaporation rate constant (atm! s7!) 100 
Keon Condensation rate constant (s7!) 100 
OHO, in Inlet water-vapor mass fraction 0.0198 
0p, in Inlet oxygen mass fraction 0.2284 
ON; in Inlet nitrogen mass fraction 0.7518 
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Table 2 

Source terms of the two-phase equations 

Equations Source term Gas diffusion Catalyst layer (CL) 

layer (GDL) 

Momentum THORR,O7 0 —jetMo,/(4F) 
THORR,H2O 0 JetMy,0/(2F) 
Mphase Eq. (8) Eq. (8) 

Energy qa 5/05. eff in, /Omett + i /0s eff 
phase Mphase x hfg Mphase x hfg 
dint hy x (Ts = Te) 0 
ORR 0 Jet x (om aa Pm) 

Charge Sj 0 jet 


respectively, represented by the following forms [14]. 
V. (—ks,eff V Ts) = —ġin + a (9) 
(pcp)u - VTe + V - (—keft V Te) = Gint + Yphase (10) 


where kg err, and kfeff, are the effective thermal conductivities 
of solid matrices and reactant fluids, respectively. In Eq. (9) 
(solid-phase temperature), the source terms ġo and int stand for 
the Joule heating and intrinsic heat transfer between the solid- 
matrices and the reactant fluids, respectively. They are provided 
in Table 2. In the fluid-phase energy equation (Eq. (10)), the 
phase-change heat transfer, phase, can be represented by the 
product of the latent heat of evaporation/condensation (hfg) and 
the interfacial mass transfer rate (phase), 1.€., 


phase = Mphase X hfg (11) 


Basically, the quantity ġphase is strongly dependent on the fluid- 
phase temperature. First, the latent heat of phase change is a 
function of the fluid temperature, i.e., Afg(Tf). Secondly, if con- 
densation occurs, the quantity itphase Will be positive; thus, 
phase > 0. That is condensation will heat up the control vol- 
ume. An increase in the fluid-phase temperature will increase 
the saturation pressure of water. Thus, the driving force for con- 
densation (pH,0 — Piso) will decrease, resulting in a decrease 
of condensation rate. A similar argument also applies to evapo- 
ration. Note that the density of the reactant fluids is expressed 
by the average value of gas-phase mixture and the liquid water, 
i.e., Of = PwS + Pg(1 — 8). 

In the catalyst layer, physically, the electrochemical reaction 
takes place at the interface of reactant fluid and catalyst under 
a fixed temperature. Therefore, the reactant fluids and the solid 
matrices have the same temperatures. 


T; = Ts (12) 
(pcp) ,U VTE + V - (—ke ett V TE) = QORR + GQ + phase (13) 


where ke eff, is the effective thermal conductivity of the catalyst 
layer [6]. The irreversiblility heating by the oxygen reduction 
reaction is expressed as GorR = Jct(¢m — $s). 


2.3. Two-phase current 


In the catalyst layer, a potential difference between the 
electrolyte and catalyst drives the transfer current that keeps 


the electrochemical reaction continuously. The current passes 
through the catalyst layer can be decomposed into two parts, 
i.e., ionic currents (im) and electronic currents (is), respectively. 
Since the electrodes are electroneutral everywhere, there is no 
charge-buildup in the catalyst layers. Thus, the charge conser- 
vation is 


V - (is + im) = 0 (14) 


These two currents interact through electrochemical reactions. 
The electrons are transferred to the electrolyte from the catalyst 
in the cathodic catalyst layer. Application of the Ohm’s law to 
the equation yields the current conservation: 


Ve (—0s eff V Ps) = jet (15) 
V(—om,etf Vm) = — jet (16) 


where the source term jet, the local transfer current density cor- 
responds to the oxygen reduction reaction in the catalyst layer. 
Os eff and Omerr, are the effective electronic and ionic conduc- 
tivities of the catalyst and electrolyte, respectively. They are 
modeled as 


Os eff = Os(1 — Ec) X Us X Dis (17) 
Om,eff = Om(1 — Ec) X Um X Pim (18) 


where vs and vm are the volume fractions of catalyst and elec- 
trolyte in the catalyst layer, respectively. pis, and pim, are the 
connectivities of catalyst and electrolyte in the catalyst layer, 
respectively. It is noted that only a long-range connection of the 
same particles stretch from the catalyst layer to the electrolyte 
ensures good conductivity. 


2.4. Multi-species transports 


The present model accounts for three species in the porous 
cathode, i.e., oxygen, water vapor and nitrogen. The species 
transports based on the Stefan—Maxwell multi-component dif- 
fusion are given by the following equations: 


A M VM 
peu: Voo, =V- p202) Di M, Vaj+ -g 
J= 
Vp 
ta-o] 0 + Sor (19) 


al M VM 
Pgu : Vono = V- Pomo) Dij M. Voj + oM 


j=l j 


Vp 
+(xj— aN + Smo (20) 


ON, = 1 — wo, — ©H,0 (21) 


The effective diffusivities of the three species in the porous cath- 
ode follow the Bruggemann model. The source terms for the 
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mass balances of oxygen and water vapor are, respectively, given 
by 


So. = MORR,O» (22) 
SHO = —Mphase (23) 


The local transfer current density is depicted by the 
Butler—Volmer correlation, i.e., 


2 
co 4a F CH0 
( = ) exp | (dm — os ) ( Z ) 
CO2 ref R Te CH) O,ref 


—4(1 — &)F 
RT; 


x exp | (dm 5 (24) 


where a is the symmetric factor and R is the universal gas con- 
stant. co, and cp,o are the concentrations of oxygen and water 
vapor, respectively. m — os is the difference of ionic potential 
and electronic potential in the porous electrode (i.e., overpoten- 
tial). It drives the transfer current density from the electrolyte 
phase to the catalyst phase that keeps the electrochemical reac- 
tion continuously. 


2.5. Boundary conditions and numerical methods 


The electrochemical and physical properties used in the cal- 
culation are given in Table 1. At the module inlet, i.e., the 
interface between the GDL and the flow channel, the con- 
centrations of gaseous species are 22.84%/1.98%/75.18% for 
O2/H20/N2, where No is considered as an inert gas and serves 
as diluents. In addition, the temperatures for cathode feeds 
are fixed at 298 K. At the interface between the GDL and the 
rib-shoulder, all species fluxes are set to zero and the temper- 
atures are fixed at 298 K. The solid potential at this interface 
is arbitrarily chosen to be zero. At the interface between the 
GDL and the CL, the electronic current and the fluxes of 
gaseous species and liquid water are continuous, while ionic 
current is set to zero. At the interface between the catalyst 
layer and the membrane, species and heat fluxes are set to 
zero and the water vapor is saturated. In addition, the total 
overpotential is used as boundary condition. For the rest of 
the boundaries, we have either insulation or symmetry condi- 
tions. 

The governing equations are numerically solved by the finite- 
element method [6—9]. Solutions are considered to be converged 
at each test condition after the ratio of residual source (includ- 
ing mass, momentum, species, and charge) to the maximum 
flux across a control surface becomes below 1.0 x 1076. A total 
of 8682 unstructured meshes were used in the computational 
domain, with a fine mesh in the reaction region. The iterations 
proceeded until the change in the calculated airflow rate between 
20 consecutive iterations was less than 0.1%. In order to achieve 
that, the necessary number of iterations varied between 300 and 
1000. A typical simulation requires about 60 min of central pro- 
cessing unit time on a Pentium IV PC (2.8 GHz, 2GB RAM) 
using Windows XP platform. 
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Fig. 2. Comparison of model results to experimental data. 


3. Results and discussion 


Before the subsequent discussion of the present results, it 
requires to validate the present model by comparing the present 
predictions with the experimental data. To this aim, a comparison 
of the cathode overpotential between the model results and the 
experimental data is made. Fig. 2 shows the experimental results 
in terms of the cathode half-cell overpotentials, and the model 
fit to the cathode side. The experimental data (solid circles) 
were determined by subtracting the activation overpotential for 
hydrogen oxidation on platinum [12] and the Ohmic loss in the 
membrane from the data obtained for the single-cell test. As seen 
in Fig. 2, the model results agree well with experimental data. 


Gas-phase velocity 
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Fig. 3. Comparison of gas-phase and liquid-phase velocity distributions. 
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Temperature Difference 
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Fig. 4. Phase temperature differences at several x stations. 


Fig. 3 shows the two-phase velocities inside the porous cath- 
ode. The gas-phase velocity vectors shown in Fig. 3(a) induce the 
gaseous mixture from the flow channel into the porous cathode. 
This is in contrast with the results obtained by the single-phase 
model [9] in which the gaseous-mixture velocity is directed from 
the porous cathode to the flow channel. As for the liquid-water 
velocity (Fig. 3(b)), itis opposed to the direction of gaseous mix- 
ture, i.e., with a direction from the cathode into the channel. As 
denoted in Fig. 3(a) and (b), the maximum velocities of gaseous 
mixture and liquid water are 7.3 x 1073 and 1.3 x 10-> ms7!, 
respectively. They occur under the elbow of the rib-shoulder 
land. 

As noted before, there are three kinds of two-phase trans- 
ports, velocity, temperature and current occurring in an acting 
porous electrode. The interactions between two phases are intrin- 
sic heat transfer, capillary force, and electrochemical reactions, 
respectively. Figs. 4-6 show the driving forces of the above 
two-phase transports. Fig. 4 shows the temperature difference 
between the solid phase and the fluid phase. In the catalyst 
layer, the reactant fluids and the solid matrices have the same 
temperature. Near the channel inlet, the solid-phase tempera- 
ture is higher than the fluid-phase temperature, while an adverse 
trend is found near the rib-shoulder. The distribution of capillary 
force along several x stations in the porous cathode is shown in 
Fig. 5. The liquid saturation is higher in the catalyst layer, due to 
water generation, than the GDL-flow channel interface. There- 
fore, the capillary force gradient form in the porous electrode 
drives the liquid water from the reaction site to the open channel. 
Fig. 6 shows the potential difference between the ionic conduc- 
tor (@;) and the electronic conductor (@,) in the catalyst layer. 
The phase potential increases along the depth (y) of the cata- 
lyst layer because of the decrease of oxygen concentration. In 
addition, the phase potential difference is higher in the region 
beneath the rib-shoulder but is lower in the region near the flow 
channel. 


Capillary Force 
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Fig. 5. Capillary forces at several x stations. 


In the porous cathode, the essential parameter for determining 
the direction of phase change is the relative humidity (RH) of 
the gaseous mixture, i.e., the ratio of partial pressure of the water 
vapor in the gaseous mixture to the saturation pressure. 

H XH20 X 
RH = eue = o Ps 
Pio Pio 


(25) 


When the relative humidity exceeds 100%, it gives rise to con- 
densation. On the other hand, evaporation occurs when the 
relative humidity is below 100% in the presence of liquid water. 
The phase-change mechanisms in a porous cathode of a PEM 
fuel cell are quiet complicated since all three parameters on 
the right-hand side of Eq. (25) vary [1]. They are competing 
in the directions of phase change and are described as follows. 
The oxygen reduction reaction in the cathodic catalyst layer will 
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Fig. 6. Overpotential (phase-potential difference) distributions at several x sta- 
tions. 
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Fig. 7. Effect of rib-shoulder temperature (Typ) on the adiabatic front of phase- 
change in the porous cathode of a PEM fuel cell: (a) Trib = 298 K; (b) Trib = 308 K; 
(c) Tib = 318 K. 


increase the molar fraction of water vapor xy,0, simply due to the 
reactant consumption. Meanwhile, it also reduces the thermo- 
dynamic pressure of the gaseous mixture pg. The combination 
of the two effects can lead to either evaporation or condensa- 
tion. In addition, the exothermic reaction in the catalyst layer 
increases the fluid-phase temperature that increases the satura- 
tion pressure Př‘ o(T®. It was found with an earlier version of 
this model that the fluid temperature can rise as high as 20 K [9], 
depending on the applied loading. This effect alone would lead 
to evaporation of liquid water. 

The net phase change is a result of the balance among 
these competing, coupled, and spatially varying mechanisms. 
Fig. 7(a)-(c) shows the phase-change characteristics inside the 
porous cathode for several rib-shoulder temperatures. In each 
plot, a phase-equilibrium front divides the porous cathode into 
two regimes, i.e., the evaporation zone and the condensation 
zone, respectively. Above the front, the liquid water evaporates 
due to higher fluid-phase temperature than the inlet tempera- 
ture. In the region below the front and above the rib-shoulder, 
the water vapor is condensed due to the decrease in fluid-phase 
temperature adjacent the cool rib-shoulder boundary. When 
the rib-shoulder temperature increases (7;jp), the condensation 


(a) 


~~Thermal-equilibrium front, 
TT; 


Thermal-equilibrium front, T,=T; 


_- Thermal-equilibrium front 


T3<T; 


Fig. 8. Effect of rib-shoulder temperature (7;;,) on the thermally adiabatic front 
in the porous cathode of a PEM fuel cell: (a) T,ip = 298 K; (b) Tin =308 K; (c) 
Trib = 318 K. 


zone reduces, as shown in Fig. 7(b) and (c). Increasing the 
rib-shoulder temperature increases the nearby fluid-phase tem- 
perature that increases the saturation pressure. Thus, it is difficult 
to condense the water vapor in this region. 

Fig. 8 shows the effect of the rib-shoulder temperature on the 
location of the thermal-equilibrium front in the porous cathode. 
On the right-hand side of the front, the fluid-phase temperature 
is higher than the solid-phase temperature. Conversely, the fluid- 
phase temperature is higher than the solid-phase temperature on 
the other site. It is seen that the thermal-equilibrium front is 
pushed leftward as the rib-shoulder temperature is increased. 


4. Summary and conclusions 


A two-phase model simultaneously considers velocity, tem- 
perature, and current in a PEM fuel cell has been developed 
in the present study. The three categories of two-phase trans- 
ports are strongly coupled in an operating fuel cell. The unique 
feature of the present model is to successfully predict the phase- 
equilibrium front and the thermal equilibrium front in a porous 
cathode of a PEM fuel cell, which allows for a more realis- 
tic spatial variation of two-phase parameters. Results show that 
increasing the rib-shoulder temperature will reduces the con- 
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densation zones. This is because the hot rib-shoulder surfaces 
increase the nearby fluid-phase temperature that increases the 
saturation pressure; thus, it is difficult to condense the water 
vapor in this region. 

In addition to an investigation of the effects of permeabil- 
ity and wetting characteristics of the porous cathode, future 
work will focus on incorporating membrane transport and fully 
coupling the anode and cathode transport, as well as on imple- 
menting more efficient numerical strategies, including parallel 
computing to allow practical use of the model for systematic 
parametric studies. 
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